{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Incremental search"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {},
   "outputs": [],
   "source": [
    "\"\"\" \n",
    "An incremental search algorithm \n",
    "\"\"\"\n",
    "import numpy as np\n",
    "\n",
    "def incremental_search(func, a, b, dx):\n",
    "    \"\"\"\n",
    "    :param func: The function to solve\n",
    "    :param a: The left boundary x-axis value\n",
    "    :param b: The right boundary x-axis value\n",
    "    :param dx: The incremental value in searching\n",
    "    :return: \n",
    "        The x-axis value of the root,\n",
    "        number of iterations used\n",
    "    \"\"\"\n",
    "    fa = func(a)\n",
    "    c = a + dx\n",
    "    fc = func(c)\n",
    "    n = 1\n",
    "    while np.sign(fa) == np.sign(fc):\n",
    "        if a >= b:\n",
    "            return a - dx, n\n",
    "\n",
    "        a = c\n",
    "        fa = fc\n",
    "        c = a + dx\n",
    "        fc = func(c)\n",
    "        n += 1\n",
    "\n",
    "    if fa == 0:\n",
    "        return a, n\n",
    "    elif fc == 0:\n",
    "        return c, n\n",
    "    else:\n",
    "        return (a + c)/2., n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Root is: 1.2414999999999783\n",
      "Iterations: 6242\n"
     ]
    }
   ],
   "source": [
    "# The keyword 'lambda' creates an anonymous function\n",
    "# with input argument x\n",
    "y = lambda x: x**3 + 2.*x**2 - 5.\n",
    "root, iterations = incremental_search (y, -5., 5., 0.001)\n",
    "print(\"Root is:\", root)\n",
    "print(\"Iterations:\", iterations)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# The bisection method"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {},
   "outputs": [],
   "source": [
    "\"\"\" \n",
    "The bisection method \n",
    "\"\"\"\n",
    "def bisection(func, a, b, tol=0.1, maxiter=10):\n",
    "    \"\"\"\n",
    "    :param func: The function to solve\n",
    "    :param a: The x-axis value where f(a)<0\n",
    "    :param b: The x-axis value where f(b)>0\n",
    "    :param tol: The precision of the solution\n",
    "    :param maxiter: Maximum number of iterations\n",
    "    :return: \n",
    "        The x-axis value of the root,\n",
    "        number of iterations used\n",
    "    \"\"\"\n",
    "    c = (a+b)*0.5  # Declare c as the midpoint ab\n",
    "    n = 1  # Start with 1 iteration\n",
    "    while n <= maxiter:\n",
    "        c = (a+b)*0.5\n",
    "        if func(c) == 0 or abs(a-b)*0.5 < tol:\n",
    "            # Root is found or is very close\n",
    "            return c, n\n",
    "\n",
    "        n += 1\n",
    "        if func(c) < 0:\n",
    "            a = c\n",
    "        else:\n",
    "            b = c\n",
    "\n",
    "    return c, n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Root is: 1.241903305053711\n",
      "Iterations: 20\n"
     ]
    }
   ],
   "source": [
    "y = lambda x: x**3 + 2.*x**2 - 5\n",
    "root, iterations = bisection(y, -5, 5, 0.00001, 100)\n",
    "print(\"Root is:\", root)\n",
    "print(\"Iterations:\", iterations)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Newton's method"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "metadata": {},
   "outputs": [],
   "source": [
    "\"\"\" \n",
    "The Newton-Raphson method \n",
    "\"\"\"\n",
    "def newton(func, df, x, tol=0.001, maxiter=100):\n",
    "    \"\"\"\n",
    "    :param func: The function to solve\n",
    "    :param df: The derivative function of f\n",
    "    :param x: Initial guess value of x\n",
    "    :param tol: The precision of the solution\n",
    "    :param maxiter: Maximum number of iterations\n",
    "    :return: \n",
    "        The x-axis value of the root,\n",
    "        number of iterations used\n",
    "    \"\"\"\n",
    "    n = 1\n",
    "    while n <= maxiter:\n",
    "        x1 = x - func(x)/df(x)\n",
    "        if abs(x1 - x) < tol: # Root is very close\n",
    "            return x1, n\n",
    "\n",
    "        x = x1\n",
    "        n += 1\n",
    "\n",
    "    return None, n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Root is: 1.241896563034502\n",
      "Iterations: 7\n"
     ]
    }
   ],
   "source": [
    "y = lambda x: x**3 + 2.*x**2 - 5.\n",
    "dy = lambda x: 3.*x**2. + 4.*x\n",
    "root, iterations = newton(y, dy, 5.0, 0.00001, 100)\n",
    "print(\"Root is:\", root)\n",
    "print(\"Iterations:\", iterations)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# The secant method"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 8,
   "metadata": {},
   "outputs": [],
   "source": [
    "\"\"\" \n",
    "The secant root-finding method \n",
    "\"\"\"\n",
    "def secant(func, a, b, tol=0.001, maxiter=100):\n",
    "    \"\"\"\n",
    "    :param func: The function to solve\n",
    "    :param a: Initial x-axis guess value\n",
    "    :param b: Initial x-axis guess value, where b>a\n",
    "    :param tol: The precision of the solution\n",
    "    :param maxiter: Maximum number of iterations\n",
    "    :return: \n",
    "        The x-axis value of the root,\n",
    "        number of iterations used\n",
    "    \"\"\"\n",
    "    n = 1\n",
    "    while n <= maxiter:\n",
    "        c = b - func(b)*((b-a)/(func(b)-func(a)))\n",
    "        if abs(c-b) < tol:\n",
    "            return c, n\n",
    "\n",
    "        a = b\n",
    "        b = c\n",
    "        n += 1\n",
    "\n",
    "    return None, n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 9,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Root is: 1.2418965622558549\n",
      "Iterations: 14\n"
     ]
    }
   ],
   "source": [
    "y = lambda x: x**3 + 2.*x**2 - 5.\n",
    "root, iterations = secant(y, -5.0, 5.0, 0.00001, 100)\n",
    "print(\"Root is:\", root)\n",
    "print(\"Iterations:\", iterations)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# SciPy implementations"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Root-finding scalar functions"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 10,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Bisection method: 1.241903305053711\n",
      "Newton's method: 1.2418965630344798\n",
      "Secant method: 1.2418965630344803\n",
      "Brent's method: 1.241896563034559\n"
     ]
    }
   ],
   "source": [
    "\"\"\"\n",
    "Documentation at\n",
    "http://docs.scipy.org/doc/scipy/reference/optimize.html\n",
    "\"\"\"\n",
    "import scipy.optimize as optimize\n",
    "\n",
    "y = lambda x: x**3 + 2.*x**2 - 5.\n",
    "dy = lambda x: 3.*x**2 + 4.*x\n",
    "\n",
    "# Call method: bisect(f, a, b[, args, xtol, rtol, maxiter, ...])\n",
    "print(\"Bisection method:\", optimize.bisect(y, -5., 5., xtol=0.00001))\n",
    "\n",
    "# Call method: newton(func, x0[, fprime, args, tol, ...])\n",
    "print(\"Newton's method:\", optimize.newton(y, 5., fprime=dy))\n",
    "# When fprime=None, then the secant method is used.\n",
    "print(\"Secant method:\", optimize.newton(y, 5.))\n",
    "\n",
    "# Call method: brentq(f, a, b[, args, xtol, rtol, maxiter, ...])\n",
    "print(\"Brent's method:\", optimize.brentq(y, -5., 5.))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## General nonlinear solvers"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 11,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "[1.24189656]\n"
     ]
    }
   ],
   "source": [
    "import scipy.optimize as optimize\n",
    "\n",
    "y = lambda x: x**3 + 2.*x**2 - 5.\n",
    "dy = lambda x: 3.*x**2 + 4.*x\n",
    "\n",
    "print(optimize.fsolve(y, 5., fprime=dy))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 12,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "    fjac: array([[-1.]])\n",
      "     fun: array([3.55271368e-15])\n",
      " message: 'The solution converged.'\n",
      "    nfev: 12\n",
      "     qtf: array([-3.73605502e-09])\n",
      "       r: array([-9.59451815])\n",
      "  status: 1\n",
      " success: True\n",
      "       x: array([1.24189656])\n"
     ]
    }
   ],
   "source": [
    "print(optimize.root(y, 5.))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 13,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "[-1.33306553]\n"
     ]
    },
    {
     "name": "stderr",
     "output_type": "stream",
     "text": [
      "/Library/Frameworks/Python.framework/Versions/3.7/lib/python3.7/site-packages/scipy/optimize/minpack.py:163: RuntimeWarning: The iteration is not making good progress, as measured by the \n",
      "  improvement from the last ten iterations.\n",
      "  warnings.warn(msg, RuntimeWarning)\n"
     ]
    }
   ],
   "source": [
    "print(optimize.fsolve(y, -5., fprime=dy))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 14,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "    fjac: array([[-1.]])\n",
      "     fun: array([-3.81481496])\n",
      " message: 'The iteration is not making good progress, as measured by the \\n  improvement from the last ten iterations.'\n",
      "    nfev: 28\n",
      "     qtf: array([3.81481521])\n",
      "       r: array([-0.00461503])\n",
      "  status: 5\n",
      " success: False\n",
      "       x: array([-1.33306551])\n"
     ]
    }
   ],
   "source": [
    "print(optimize.root(y, -5.))"
   ]
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 3",
   "language": "python",
   "name": "python3"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 3
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython3",
   "version": "3.7.0"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 2
}
